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Spray  Formation:  Three  Dimensional  Liquid 
Jet  Breakup  Due  to  Surface  Tension 

In  spray  combustion,  a  liquid  fuel  is  atomized  and  distributed  in  the 
combustion  chamber  to  enhance  the  combustion  process.  Atomization  is  achieved 
by  forcing  a  jet  of  liquid  fuel  through  a  fuel  injector  and  into  a  gas  whereby  it 
is  broken  up  into  a  large  number  of  minute  drops.  The  atomization  process  is 
critical  to  the  combustor’s  performance.  The  focus  of  this  study  is  the  initial 
region  of  the  liquid  fuel  jet,  as  it  begins  to  break  up  under  the  influence  of 
fluid  dynamic  and  surface  tension  forces.  Liquid  jet  breakup  is  difficult  to  study 
experimentally  due  to  the  small  spatial  and  temporal  scales  involved.  However, 
the  flow  is  governed  by  a  well-defined  set  of  equations,  and  thus,  it  can  be 
studied  by  mathematical  means. 

The  present  work  involves  analysis  and  numerical  simulation  of  the  initial 
stage  of  liquid  jet  breakup  which  is  a  precursor  to  atomization.  The  objectives 
of  the  present  work  are  to  improve  the  understanding  of  the  instability 
mechanisms  in  the  initial  region  of  the  liquid  jet  and  to  explore  possible 
methods  for  exploiting  surface  instabilities  to  improve  the  atomization 
characteristics  of  the  fuel  injection  process.  In  this  first  year  of  a  two  year 
study,  the  numerical  methods  needed  to  investigate  the  liquid  jet  breakup  have 
been  evaluated  and  assembled.  Initial  calculations  of  two  phase  flows  have  been 
made.  Also,  mathematical  analysis  has  been  performed  on  one  aspect  of  jet 
breakup  and  is  in  progress  on  another  aspect.  These  topics  are  discussed  below. 

Numerical  simulation  of  the  liquid  jet  is  performed  by  integrating 
numerically  the  discretized  Navier-Stokes  equations  for  an  incompressible,  variable 
density  fluid.  The  simulation  is  done  on  a  fixed  grid  in  two  or  three  space 
dimensions  and  time.  As  originally  planned,  the  numerical  algorithm  used  to 
integrate  the  equations  is  a  finite  volume,  fractional  steps  method.  However,  the 
discrete  density  jump  imposes  restrictions  on  specific  details  of  the  method,  and 
a  significant  amount  of  time  was  invested  in  a  study  of  Navier-Stokes  integration 
methods. 

For  a  fundamental  study  of  liquid  jet  dynamics,  it  was  felt  that  the 
algorithm  must  have  good  accuracy  and  must  guarantee  mass  and  momentum 
conservation  at  the  sharp  density  discontinuity  at  the  liquid/gas  interface. 
Aerodynamics  research  has  shown  that  nonconservative  metnods  will  give 
incorrect  motion  of  discontinuities,  which  would  be  intolerable  in  the  present 
work.  This  dictates  that  the  true  momenta,  pUj,  rather  than  the  velocities,  Uj, 
be  conserved  variables.  For  flows  with  large  density  differences  between  liquid 
and  gas  phases,  it  then  becomes  essential  that  the  density  implicitly  associated 
with  the  momenta  must  be  consistent  with  the  density  computed  by  the  mass 
conservation  equations,  since  small  inconsistencies  could  produce  large  fluctuations 
in  computed  velocities.  The  only  discretization  scheme  which  could  be  devised 
that  would  satisfy  the  conservation  and  consistency  requirements  was  a  cell 
centered  scheme.  The  widely  used  staggered  grid  discretisation  scheme  uses 
different  cells  within  the  grid  for  conserving  mass  and  the  momenta,  and  as  a 
result,  is  essentially  incompatible  with  the  consistency  requirement.  A  decision 
to  develop  a  new  algorithm  was  made  after  a  thorough  assessment  of  existing 
Eulerian  methods  revealed  that  none  of  the  methods  appropriate  for  the  work  (as 
described  in  the  literature)  satisfied  momentum  conservation  in  the  presence  of 
density  discontinuities. 
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Thus,  a  method  which  uses  a  nonstaggered  grid  was  developed.  The  density 
and  momenta  are  stored  at  cell  centers,  and  the  pressure  is  stored  at  cell 
corners.  Otherwise,  it  is  similar  to  the  methods  envisioned  in  the  original 
proposal.  The  present  method  is  second  order  in  time  and  second  order  in 
space  with  optional  fourth  order  accuracy  on  convective  terms.  It  uses  the 
fractional  steps  concept  to  satisfy  the  divergence-free  constraint  while  maintaining 
good  efficiency. 

An  interface  tracking  method  is  used  in  conjunction  with  the  conservation 
equations  to  insure  that  the  density  discontinuity  is  resolved  sharply.  This 
approach  is  second  order  accurate,  which  is  felt  to  be  the  minimum  required  for 
studying  stability.  The  Simple  Line  Interface  Calculation  (SLIC)  method 
originally  envisioned  for  this  work  is  first  order,  and  computed  results  reveal 
certain  characteristics  (e.g.,  the  truncation  of  small  waves  on  an  interface  nearly 
parallel  with  grid  lines)  that  would  hamper  a  study  of  interface  instability 
mechanisms.  Surface  tension  is  included  via  the  interface  tracking  method  and 
the  momentum  equations.  In  general  the  surface  tension  method  gives  good 
results.  However,  there  was  a  problem  in  the  surface  tension  procedure  for  some 
shallow-angle  waves  which  caused  the  method  to  become  unstable.  The  problem 
was  traced  to  an  approximation  in  the  fractional  steps  procedure,  and  it  has 
been  partially  corrected.  Some  further  work  is  required,  but  this  is  believed  to 
be  straightforward.  The  results  given  below  avoid  this  problem. 

A  computer  code  which  employs  the  above  described  methods  was  written 
and  is  being  evaluated  on  some  sample  problems.  Single  phase  flow  in  a  cavity 
with  a  driven  lid  was  computed,  as  this  is  a  ‘standard’  test  case  for  Navier- 
Stokes  methods.  Accurate  results  at  several  Reynolds  numbers  have  been 
computed.  Two-  and  three-dimensional  cavity  flow  calculations  have  been 
performed  successfully.  The  remaining  calculations  are  of  two-dimensional 
axisymmetric  flows. 

In  a  uniform  flow  without  gravity  or  surface  tension,  a  surface  wave  should 
convect  indefinitely  without  distortion.  The  ability  to  simulate  this  distortion- 
free  convection  is  critical  to  the  ability  to  predict  accurately  the  growth  rates 
of  waves  on  the  liquid  surface.  The  density  jump  is  used  to  identify  the 
interface.  (The  density  contours  show  a  smeared  interface,  but  this  is  an  artifice 
of  the  plotting  package.)  Figure  1(a)  shows  the  initial  condition  of  the  wave. 
Figures  l(b,c)  show  the  wave  after  it  has  traveled  one  half  and  two  wavelengths, 
respectively,  and  the  fidelity  of  convection  prediction  can  be  seen.  Although  this 
test  case  is  physically  uninteresting,  it  is  a  challenging  numerical  test  because 
there  is  a  10  to  1  jump  in  density  and  momentum  across  the  interface  which 
must  be  convected  through  all  cells  between  the  trough  and  crest  of  the  wave. 
The  present  method  gives  good  results  for  this  flow,  and  it  shows  none  of  the 
‘truncation’  problems  associated  with  the  SLIC  method. 

The  Rayleigh  instability,  which  causes  a  liquid  cylinder  with  surface  tension 
to  break  up  into  drops,  has  also  been  simulated.  The  liquid  cylinder  is  initially 
at  rest,  but  its  surface  is  given  a  small  sinusoidal  displacement.  Surface 
tension  induces  a  pressure  disturbance  which  causes  long  waves  to  be  amplified 
and  short  ones  to  decay.  As  above,  there  is  a  10  to  1  jump  in  density  across 
the  interface.  Figures  2(a,b)  show  the  density  contours  and  velocity  vectors  for 
a  long  wave,  and  the  velocities  indicate  that  the  disturbance  is  growing.  The 
decay  of  a  short  wave  is  shown  in  Figures  3(a,b).  Quantitative  evaluation  of 
the  method’s  accuracy  for  this  problem  has  not  yet  been  performed,  but  it  is 
planned  for  the  near  future. 


The  above  test  cases  have  done  nothing  to  increase  our  understanding  of  the 
physics  of  the  problem.  However,  they  do  show;  that  the  numerical  method  can 
treat  the  sharp  density  jump  and  surface  tension  at  the  interface,  even  for  short 
wavelengths.  Since  the  object  of  the  work  is  to  investigate  the  physics,  it  is  felt 
that  the  use  of  an  accurate  method  is  essential.  These  test  cases  indicate  that 
the  desired  accuracy  has  been  achieved. 

At  present,  the  work  is  directed  at  completing  the  evaluation  of  the 
computer  code  and  making  any  needed  modification.  The  ability  to  simulate  the 
conical  liquid  jet  should  be  achieved  early  in  the  second  year,  at  which  time  the 
work  will  focus  on  the  instability  mechanisms  of  fuel  sprays. 

Other  aspects  of  interface  flows  have  also  been  investigated.  These  are 
being  pursued  at  a  low  level  of  effort  and  are  intended  to  complement  the  main 
computational  thrust  of  the  work.  Inviscid  stability  analysis  of  the  cylindrical 
liquid  sheet  (of  constant  mean  radius)  has  been  used  by  previous  researchers  to 
explain  some  conical  fuel  spray  phenomena.  In  the  present  work  this  analysis 
was  extended  to  a  liquid  sheet  with  a  radius  which  increases  in  time.  The 
analysis  should  more  closely  approximate  the  physics  of  conical  sprays  than  the 
constant  radius  case.  However,  it  was  found  that  the  stability  mechanisms 
which  are  present  in  the  constant  radius  case  are  unaffected  by  the  rate  of 
change  of  radius.  Inviscid  Lagrangian  computational  methods  were  also 
investigated.  These  methods  are  well  suited  to  two-dimensional  simulations  and 
are  being  considered  as  a  means  of  checking  the  Navier-Stokes  simulations.  An 
Orr-Sommerfeld  analysis  of  a  shear  layer  with  variable  density  and  viscosity  has 
been  started  but  not  completed.  If  successful,  this  work  will  provide  guidance 
on  viscous  instability  mechanisms,  and  it  may  help  explain  Reynolds  number 
effects  on  atomization. 

Progress  this  year  has  been  slightly  slower  than  planned  due  to  difficulties 
with  the  Navier-Stokes  integration  algorithm.  The  fundamental  shortcomings  of 
previous  methods  for  flows  with  discontinuous  density  were  not  fully  anticipated. 
It  is  believed  that  the  time  spent  to  resolve  the  difficulties  which  were  uncovered 
will  yield  more  accurate  results  for  spray  flows.  This  delay  should  not  prevent 
the  research  from  achieving  its  stated  goals. 

Drs.  Robert  E.  Childs  and  Nagi  N.  Mansour  have  worked  on  this  project. 
No  degrees  have  been  awarded  in  conjunction  with  this  work.  No  publications 
have  yet  resulted  from  this  work,  but  it  is  planned  that  an  abstract  will  be 
submitted  to  the  26th  ALA  A  Aerospace  Sciences  Meeting,  to  be  held  in  January 
1988.  Mr.  Dan  Bulzan  of  the  NASA  Lewis  Research  Center  has  been  contacted 
regarding  this  work.  General  objectives  were  covered.  No  new  discoveries  or 
patents  have  resulted  from  this  work. 
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Figure  1.  Density  contours  show  essentially  distortion-free  convection  of 
a  simple  sinusiodal  wave  in  a  uniform  flow  field. 


Figure  2.  Density  contours  and  velocity  vectors  for  a  growing  wave  on 
a  liquid  cylinder  with  surface  tension. 
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Figure  3.  Density  contours  and  velocity  vectors  for  &  decaying  wave  on 
a  liquid  cylinder  with  surface  tension. 
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